Dynamical transition, hydrophobic interface, and the temperature dependence of 

electrostatic fluctuations in proteins 
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Molecular dynamics simulations have revealed a dramatic increase, with increasing temperature, 
of the amplitude of electrostatic fluctuations caused by water at the active site of metalloprotein 
plastocyanin. The increased breadth of electrostatic fluctuations, expressed in terms of the reor- 
ganization energy of changing the redox state of the protein, is related to the formation of the 
hydrophobic protein/water interface allowing large-amplitude collective fluctuations of the water 
density in the protein's first solvation shell. On the top of the monotonic increase of the reorganiza- 
tion energy with increasing temperature, we have observed a spike at 220 K also accompanied by a 
significant slowing of the exponential collective Stokes shift dynamics. In contrast to the local den- 
sity fiuctuations of the hydration-shell waters, these spikes might be related to the global property 
of the water solvent crossing the Widom line. 

PACS numbers: 87.14.E-, 87.15. N-, 87.15. Pc 
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I. INTRODUCTION 

Dynamical transition has been observed in many hy- 
drated biopolymers, including proteins, DNA, and RNA 
[U, 0, H, 3 ■ In amounts to a sharp change in the temper- 
ature slope of mean-squared atomic displacements of the 
biopolymer atoms at the temperature usually observed in 
the range Ttr = 200 — 230 K. While the microscopic origin 
of this dynamical transition is still debated Pi fi S B 9 ' 
an important open question is how the existence of this 
universal property of hydrated biopoly mers [lO| affects 
their physiological activity [l|, 0, llO, llll ■ 

Electrostatics is significant to the catalytic action of 
enzymes [i2j- Therefore, a link between protein's dy- 
namical transition and enzymatic activity may exist in 
some property characterizing electrostatics at the active 
site. It is currently well established that the dynamical 
transition is not observed in dry proteins, and its exis- 
tence is universally attributed to the interaction of water 
with the protein interface. A property sensitive to the dy- 
namical transition needs to connect water's electrostatics 
to protein's active site. Here, we consider one such pa- 
rameter which critically affects barriers of protein redox 
reactions, the reorganization energy of electron transfer 

The reorganization energy A of electronic transitions 
between proteins characterizes the breadth of thermal 
fluctuations of the energy gap AE between the donor 
and acceptor energy levels 
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Here, SAE is the fluctuation of the energy gap AE and 
the inverse temperature /3 = l/{kBT) corrects for the 
proportionality of the variance to temperature following 
from the fluctuation-dissipation theorem [l4|. The reor- 
ganization energy A is then typically a weak function of 
temperature when measured for electronic transitions in 
molecular polar solvents [ISj . 

Experimentally accessible reorganization energy of in- 



terprotein electron transfer [l6| characterizes the cou- 
pling of the energy levels of both the donor and acceptor 
to the thermal bath. For long-distance electron transfer, 
most common in biological energy chains, A can be split 
into a sum of individual, donor and acceptor, components 
and a Coulomb correction. Since these individual com- 
ponents mostly characterize the physics of the problem, 
our focus here is on the electrostatic fluctuations at the 
active site of a single protein. 

Electron transfer changes the redox state of the protein 
and thus the partial atomic charges of the active site. The 
electrostatic interactions of these charge differences with 
the potential of the hydrating water (/>„, j- at atomic sites 
j contribute to the Coulomb shift AE^ which is a part 
of the overall donor-acceptor energy gap: 
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Here, the sum runs over the atoms of the active site. The 
variance of this Coulomb energy gap calculated for the 
charges Aqj of the active site of a single protein is what 
is studied in this paper. The water reorganization energy 
is then defined as 
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The dynamical dimension of the problem is character- 
ized by the normalized Stokes shift correlation function 

m 

S^lt) = {SAE2{t)6AE^m/{i6AE^iO)f), (4) 

where angular brackets denote an ensemble average. The 
common form of Sm{t) is dense polar liquids includes a 
fast one-particle component with a Gaussian decay fol- 
lowed by exponential (or stretched exponential) decay 
describing collective solvent dynamics [131 
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FIG. 1: Stokes shift correlation function of PC (Ox) at dif- 
ferent temperatures indicated in the plot 



Here, tg and te are, respectively, the Gaussian and ex- 
ponential relaxation times and Aq quantifies the relative 
weight of single-particle dynamics in the reorganization 
energy; Pe is the stretching exponent. 

The main result of the Molecular Dynamics (MD) 
simulations presented here is to show that the reorga- 
nization energy \^{T) rises significantly with tempera- 
ture to a value much exceeding both the common es- 
timates of this parameter for reactions involving small 
redox molecules, [13 and previous estimates for protein 
electron transfer [l^l. We associate this increase with 
the formation of the hydrophobic interface allowing large- 
amplitude fluctuations of the local water density. We also 
show that both the long-time exponential relaxation time 
[te in Eq. ([5])] and the collective part of the reorganiza- 
tion energy Xe = (1 — Aq)Xw pass through peaks at the 
temperature of dynamical transition Ttr = 220 K. This 
special temperature (at atmospheric pressure) has been 
previously associated with a thermodynamic singularity 
in the phase diagram of bulk water [18|, |19| . It has also 
been recently suggested that a transition from fragile to 
strong dynamics of hydrated biopolymers occurs at the 
same temperature 0, [20, [21| . 



II. MD SIMULATIONS 

MD simulations reported here have been done for the 
redox metalloprotein plastocyanin (PC) from spinach ac- 
cording to the simulation protocol described in our pre- 
vious publication [22|. PC is a single polypeptide chain 
of 99 residues forming a /3-sandwich, with a single copper 
ion ligated by cysteine, methionine, and two histidincs. 
The protein's active site in our analysis is composed of 
a copper ion and four atoms (two nitrogens and two sul- 
furs) coordinating it. The partial charges on these atoms 
in both reduced (Red) and oxidized (Ox) states can be 
found in Ref. HI. 

The initial configuration of PC was taken from a pro- 
tonated X-ray crystal structure with a 1.7 A resolution 
(PDB: lag6, [2^). First, the initial protein configura- 
tion was minimized in vacuum using the conjugate gra- 
dient method for 10^ steps to remove any bad contacts. 
Then, the system was solvated in an octahedral box with 



N^ = 5886 TIP3P molecules [25|, providing at least two 
solvation shells around the protein. The protein was sim- 
ulated in the Ox state with a total charge of —8 and in 
the Red state with the total charge of —9. In both cases, 
eight or nine sodium ions were added to neutralize the 
system, as is required for the Ewald summation. After 
adding the water and counterions, the system's energy 
was minimized for another 10^ steps while the protein 
was allowed to relax and the water and sodium atoms 
were positionally constrained. Finally, the entire system 
was additionally minimized for 10^ steps. 

Following minimization, the system was heated in a 
NVT ensemble for 30 ps from K to the desired tem- 
perature. Temperature equilibration was followed by a 
2 ns density equilibration in a NPT ensemble at P = 1 
atm. This equilibrated structure was then used for 20 
individual simulations of the Ox state and 7 simulations 
of the Red state of PC to create 10 ns long trajecto- 
ries. Temperatures T was varied from 100 to 300 K 
at constant volume V and constant number of water 
molecules N^ ■ The total simulation time was 324 ns and 
required 6.9 CPU years, while only 270 ns were used for 
the production data analysis which lasted another 2.2 
CPU years. The timestep for all MD simulations was 
2 fs, and SHAKE was used to constrain bonds to hy- 
drogen atoms. Constant temperature and pressure sim- 
ulations employed Berendsen thermostat and barostat, 
respectively [2Q] . The long-range electrostatics were cal- 
culated using a smooth particle mesh Ewald summation 
with a 9 A limit in the direct space sum. 



III. RESULTS 

Most MD results reported here have been obtained 
from configurations in equilibrium with the Ox state of 
PC; the Stokes shift data were collected from both Ox 
and Red equilibrium trajectories as discussed below. The 
reorganization energies Xw{T) of PC(Ox) state were ob- 
tained from MD trajectories at different temperatures T 
and fixed observation window Tobs = 1 ns. More specif- 
ically, the reorganization energy is calculated from the 
variance of the Coulomb energy gap [Eq. ([1])] by sliding 
a 1 ns observation window along a longer MD trajectory 
and averaging over the results of the variance calcula- 
tions on each window. The average (A£'^)obs required 
to calculate the variance is not a global average but is 
obtained separately from each observation window. This 
approach to the calculation of averages is analogous to 
a laboratory procedure with a fixed resolution and is re- 
quired for studies of systems with broad distributions of 
relaxation times ^^. In case of proteins, a subset of 
nuclear motions is always frozen on the simulation time- 
scale and so both specifying the observation window and 
keeping it constant for all measurements is significant in 
maintaining consistent conditions for collecting the data. 

Fits of the simulated Stokes shift functions to Eq. ([5]) 
are shown in Fig. [TJ Two features are most prominent 
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FIG. 2; Panel (a): Points are the water reorganization ener- 
gies Au, from MD simulations of PC (Ox) and the solid line 
shows the fit of the simulation data to Eq. (O with the fit- 
ting parameters: Xq = 0.39 eV, Acq = 0.87 -I- 0.0084 * T eV, 
Toba/To — 7350, and Ea = 1867 K; the dashed line assumes 
temperature- independent Acq. The inset shows the exponen- 
tial part of the reorganization energy related to collective wa- 
ter fiuctuations. The closed diamonds refer to half of the 
Stokes shift [Eq. (|8])] and the open diamonds show the linear- 
response reorganization energy A"" obtained from Eq. ((9}. 
Panel (b): The points refer to the exponential relaxation time 
of the Stokes shift correlation function obtained from MD of 
PC (Ox) and the inset shows the Gaussian amplitude Ac in 
Eq. ©. 



there: the increase of the relative importance of the 
collective solvent dynamics with increasing temperature 
(decrease of Aq in Eq. ^\ and the appearance of a peak 
in the exponential relaxation time at Ttr = 220 K (Fig. 
[2}d). The exponential part of the reorganization energy 
\e also shows a peak at the same temperature (inset in 
Fig.EK). 

Overall Au,(T) strongly increases from a value typi- 
cal for short MD simulations of proteins [12] to a much 
larger value at higher temperatures (Fig. [2^ and Table 
U]). The temperature of the onset of the Au,(r) rise is 
much below Ttr, at about 150 K commonly associated 
with the onset of rotation of methyl groups of protein's 
side chains [a, l27|. This onset temperature is however 
depends on the observation window. Since the relax- 
ation times of the protein are widely different, the rise of 
\w {T) is caused by the appearance of a particular relax- 
ation mode in the observation window, methyl rotations 
in this case. However, we believe that the underlying 
picture is more complex and the main rise of Au,(T) is 
caused not by methyl rotations, but by a more collec- 
tive mode coupled to the solvent interfacial translations 
|28ll29J (see below). In fact recent extensive simulations 
of the mean-squared atomic displacements of myoglobin 
[27| have reveled two breaks in the temperature slope: 
the first break at 150 K related to methyl (anharmonic) 
rotations followed by a stronger solvent-induced break at 



FIG. 3: Average number of water molecules in PC's first sol- 
vation shell (a) and its variance (b) vs temperature. The 
dashed line in (b) shows a linear regression through the points 
and the dotted line in (a) connects the simulation points. 



TABLE I: Properties of hydrated plastocyanin (Ox) from 
MD simulations. 
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"Water reorganization energies (in cV), obtained with r^hs = 1 ns 
observation window. 

''Diffusion coefScients of TIP3P water averaged over all molecules 
in the simulation box (in A^/ps). 

■^Exponential relaxation time (in ps) of the Stokes shift correlation 
function in Eq. JSj. 



220 K. 

The appearance of a relaxation mode in the observa- 
tion window restores the statistical ergodicity for that 
particular mode. The non-ergodic rise of XwiT) to its 
equilibrium value Aca(T'), also seen for model charge- 
transfer chromophores [30|, can be described by impos- 



ing a step-wise frequency filter on the spectrum of Stokes 
sliift fluctuations [iq 



A^(r) = 2Aeq(T) 
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Here, Sw{^) is the Fourier transform of the Stokes shift 
correlation function in Eqs. ^ and ([5]). In order to pro- 
vide a physically transparent form for A^(T) one can 
consider an effective single-exponential Debye relaxation, 
instead of several relaxation modes, to characterize col- 
lective nuclear motions coupled to the Stokes shift dy- 
namics. This procedure leads to the following simple 
relation 

A^(r) = AG + (Aeq(T)-AG)(2/7r)arctg[Tobs/T(r)] , (7) 

where the effective Debye relaxation time is given by the 
Arrhenius law, t{T) = roexp[/3_Ba]- 

The Gaussian component of the solvent reorganization 
energy, related to ballistic water motions [13| , is normally 
reasonably temperature- independent 3011. On the other 
hand, the temperature decrease of Ac{T) ~ Xc/XcqiT) 
in the fit of the Stokes shift function (inset in Fig. [^b) 
clearly points to the equilibrium reorganization energy 
increasing with temperature. From the anticipated rela- 
tion of Acq(T) with the variance of the number of parti- 
cles in the first solvation shell, which linearly grows with 
temperature (see below), we have attempted a linear tem- 
perature dependence of Aeq(T) to fit the MD data to Eq. 
([7]) . The result is shown by the solid line in Fig. [2^, and 
it is not much different from the fit using a temperature- 
independent Acq (dashed line in Fig. [2^) • We also note 
that since our simulation length obviously cuts some slow 
nuclear modes off, we have not used Ag{T) from the 
fits of the Stokes shift correlation functions to calculate 

Acq(T). 

The activation energy Ea of the effective Debye mode 
obtained from the fit, Ea — 1867 K, points to a sec- 
ondary /3-relaxation mode creating fluctuations of the 
electrostatic potential, in contrast to the primary a- 
relaxation of the water/protein system with a commonly 
much higher activation barrier [3l|, [S^] . This activation 
energy is also lower than (3 relaxation of aqueous mix- 
tures with the activation energy of the order of 5.5 x 10'^ 
Ki. 

Also shown in Fig. [2^ (closed diamonds) is the reorga- 
nization energy from the Stokes shift obtained from the 
difference of average Coulomb energy gaps in Ox and Red 
states: 
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For water fluctuations following linear response one ex- 
pects the reorganization energy from the variance [Eq. 
(|3|)] to be connected to the reorganization energy from 
the two flrst moments [Eq. ([5])] by the following relation: 



In Eq. ^ we have stressed that the averages are taken 
over the configurations in equilibrium with PC (Ox) and 
AEp is the Coulomb interaction energy of the difference 
charges of the active site [Ag^, Eq. ^] with the remain- 
ing partial charges of the protein matrix. 

When a rigid molecule is solvated and the intramolecu- 
lar energy gap does not fluctuate, the second correlator in 
Eq. ^ is zero. One arrives then at the standard expec- 
tation of the linear solvation theories that two routes to 
the reorganization energy, from the second cumulant [Eq. 
dS])] and from two first cumulants [Eq. ([5])], are equiva- 
lent [13, [la| . Since the protein matrix fluctuates itself, 
the cross-correlation in principle needs be taken into ac- 
count, and it turns out to be negative [33|. However, 
when cross-correlation term is subtracted from A^* in Eq. 
Q (open diamonds in Fig. [2^), the result is still signifl- 
cantly below the reorganization energy from the variance 
[Eq. ([3])]. We therefore observe here a severe breakdown 
of linear solvation. 

What Fig. [2^ in fact indicates is that the two deflni- 
tions of the reorganization energy converge at low tem- 
peratures with the reorganization energy from the vari- 
ance deviating signiflcantly upward above T ~ 200 K. 
This observation implies that fast water's modes cou- 
pled to electrostatic fluctuations, presumably librations, 
which are still unfrozen at low temperatures, follow the 
expectations of the linear response theories. On the con- 
trary, a slower collective mode, which appears in the ob- 
servation window at higher temperatures and gives rise 
to the gigantic reorganization energy, does not follow the 
linear response. The cross-correlation does not restore 
the linear response, in contrast to an earlier observation 
made for a water-exposed tryptophan residue [33||. The 
low value of the cross-correlation physically implies that 
the elasticities of the protein and water are drastically 
different and their electrostatic fluctuations are mostly 
decoupled. From that perspective, this correlation de- 
coupling should hold for any solute/solvent combination 
with a signiflcantly different rigidity. 

The fluctuations of water's electrostatic potential at 
the active site can generally be traced back to two weakly 
correlated nuclear modes in polar liquids, the orienta- 
tional polarization and the local density [34| . In order 
to clarify the origin of the dramatic rise of the reorga- 
nization energy, we have looked at two additional corre- 
lation functions characterizing the density and orienta- 
tional manifolds of the water molecules in the protein's 
first solvation shell. A water molecule is defined as to 
belong to the first solvation shell if its oxygen atom is 
within 2.87 A distance from the protein van der Waals 
surface. 

The density manifold is characterized by the fluctua- 
tion of the number of particles Ni{t) in the flrst solvation 
shell, 
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Further, the orientational manifold is described by the 
fluctuations of the total dipolc moment M(i) of the water 
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FIG. 4; {{SM)'^)/Ni{T) vs temperature. The inset shows the 
second-rank orientational order parameter P2{T) [Eq. (|12|l ]. 



dipoles in the first solvation shell, 

CM(i) = (A^/)"'(^M(t).5M(0)), 



(11) 



where Nj{T) — (Nj) is the average number of waters in 
the first solvation shell. In CN{t) and CM{t), the fluctu- 
ations 6Nj{t) and (5M(t) denote the deviations from the 
corresponding average values. The variances were calcu- 
lated on the 1 ns observation window by using the same 
procedure as for the reorganization energy calculations. 
The temperature dependences of the average and vari- 
ance of the number of waters in the first solvation shell 
(Fig. [2]) are indicative of the formation of the hydropho- 
bic protein/water interface with increasing temperature. 
The average Nj (T) is generally a decaying function (Fig. 
|3^) , and the slope of this decay becomes sharper above 
the transition temperature Ttr (see below). The decrease 
in the density of water at the interface allows stronger 
density fluctuations (Fig. [SJd) and it is this regime of large 
interfacial density fluctuations that is a signature of hy- 
drophobic solvation [3^. In this regime, one-particle ex- 
changes of water molecules between the surface and the 
bulk [2!|| combine into large-scale collective density waves 
producing signiflcant modulations of the electrostatic po- 
tential reflected in Au,(T). This thermal noise of hy- 
drophobic surfaces is also reflected in a well-documented 
increase of protein's heat capacity upon unfolding, in- 
dicative of an increased breadth of the energy fluctua- 
tions Hi, [33. 

The interfacial density fluctuations originate from the 
exchange of waters between the hydration shell and the 
bulk. These fluctuations can be represented as bind- 
ing/unbinding events at the protein surface [38| with the 
resulting equilibrium reorganization energy Aeq(T) scal- 
ing linearly with the variance of the number of particles in 
the hydration shell: A^(r) = a + b{{SNiy{T)), where co- 
efficients a and b are weak functions of temperature. This 
expectation, used in the solid-line fit in Fig. [2^, is corrob- 
orated quite well given the linear scaling of {{6Ni)'^{T)) 
with temperature (Fig. [5]d). A fairly significant temper- 
ature rise of Acq(T) (see the fitting parameters in Fig. 
[2]) also indicates a substantial density component in the 
overall reorganization energy at ambient conditions, in 
contrast to a 20-30% contribution for small solutes in 
dense polar solvents [SJ]. We therefore conclude that 
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FIG. 5: Exponential relaxation time t^' extracted by fitting 
the correlation function CM{t) from Eq. 1)11^ to Eq. ([Sjl. The 
points are the simulation results in three ranges of tempera- 
ture where they are fitted to the Vogel-Fulcher temperature 
law (at highest temperatures) and to Arrhenius laws (at in- 
termediate and lowest temperature ranges). The solid lines 
show the results of the fits. 



the contribution to Aoq(T) from density fluctuations is 
signiflcantly magnifled by the soft and flexible nature of 
the hydrophobic protein/water interface resulting in a gi- 
gantic magnitude of the overall reorganization energy far 
exceeding half of the Stokes shift [Eqs. ^ and ^]. 

The orientational fluctuations of the first-shell dipoles 
do not show a resolvable correlation with the reorganiza- 
tion energy (Fig. |3|) . The variance of the first-shell dipole 
moment grows with rising temperature, in accord with 
a general expectation of increased softness of the solva- 
tion shell, but does not show an obvious correlation with 
Xw{T). There is a weak maximum at Ttr for (((5M)^), but 
it is hard to assess from our data whether this is another 
refiection of the same spike seen for the collective part of 
the reorganization energy A^; in Fig. [2^. 

The inset in Fig.|3]shows the second-rank orientational 
order parameter 



P2 




(12) 



Here, e^ and f j are the unit vectors of the dipole mo- 
ment and position of molecule j which belongs to the first 
solvation shell, and P2{x) is the second Legendre polyno- 
mial. The low-temperature portion oip2{T) is practically 
constant showing a slight preferential orientation of the 
water molecules parallel to the interface. This type of or- 
dering has been previously observed at interfaces of non- 
polar substances and proteins with water [33, Sol ■ This 
preferential ordering decays with increasing temperature 
resulting in essentially random, on average, orientations 
of water dipoles in the hydration shell. The fairly large 
amplitude of the dipole moment fluctuations is therefore 
most likely caused by the density fluctuations. 

This assessment is supported by the data for exponen- 
tial relaxation times of C^it) and CA/(t) obtained by 
fitting these correlation functions to Eq. ([5]). When both 
exponential relaxation times are fitted to Arrhenius laws, 
they produce activation energies of 1389 K and 2076 K, 



respectively, in a close range with the activation energy 
of 1867 K obtained from the fit of A^(T) to Eq. ^. We 
note that this activation barrier is consistent with the 
activation enthalpies of 1400-2400 K obtained by a vari- 
ety of techniques for /3-fluctuations of hydrated proteins 
[32J which are considered to be slaved by /?- fluctuations 
of the hydration shell [32, |4i|. One also needs to keep 
in mind that an average Arrhenius slope actually hides a 
fairly complex behavior. Figure [5] shows the exponential 
relaxation time of Cmit) vs inverse temperature. The 
low-temperature portion of the data (triangles) is well 
approximated by a non- Arrhenius Vogel-Fulcher temper- 
ature law. This is followed by what can be characterized 
as a fragile-to-strong crossover followed by yet another 
break in the Arrhenius slope at ~ 160 K. This picture is 
consistent with two breaks in the slope seen in the simula- 
tions of mean-squared atomic displacements of myoglobin 
[27| . where the lowest-temperature break was associated 
with the onset of methyl group rotations. The results for 
exponential relaxation times of Cn (t) are more scattered 
and we could not reach an equally informative conclusion 
except for the average Arrhenius slope. 



IV. DISCUSSION 

Many alternative explanations have been sought for 
the observed dynamical transition in biopolymers \w^ . 
Given that the transition is not observed for dry protein 
samples, the possible scenarios are limited to either the 
protein/water interface or to a bulk property of water. 
The recent observations, from neutron scattering mea- 
surements, of the fragile-to-strong crossover in the d 
namics of partially hydrated protein powder samples 
point to the second (bulk water) scenario. The crossover, 
also seen in the recent simulations [g, [2l|, can be con- 
nected to the bulk water crossing the Widom line, i.e. the 
line of maximum cooperativity of the water fluctuations 
|18| . On the other hand, other recent experimental data 
on quasielsatic neutron scattering, dielectric relaxation 
[8| , and conductivity [7,] of hydrated proteins have not re- 
vealed any special points in the corresponding relaxation 
times around the temperature of dynamical transition. 
These latter data report the temperature dependence of 
the primary a- relaxation of the protein/ water system and 
therefore these authors have concluded that the observed 
dynamical crossover [Sj should be attributed the appear- 
ance of a secondary relaxation in the observation window 
atT>rt, 1,113. 

In addition, recent observations of the dynamical tran- 
sition in DNA and RNA [13, El have clearly shown that 
this property is not unique to a peptide-based polymer. 
These findings again re-emphasize the notion that either 
a bulk property of water or some generic property of the 
interface, not much sensitive to the details of the macro- 
molecular structure, are responsible for the transition. 
Our data in fact suggest that both bulk and interfacial 
views need to be invoked to explain different facets of the 
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FIG. 6: Diffusion coefRcient (open circles) from the present 
simulations and atomic mean-squared displacements of myo- 
globin measured experimentally [32| (small up-triangles) and 
obtained from MD simulations [27| (small squares). Closed 
diamonds show the change in the number of particles in the 
first solvation shell (Fig. |3]). All parameters are normalized 
to their corresponding values at 300 K. 



problem, but the interface aspect has a dominant effect. 

We have shown that the dramatic rise of the reorga- 
nization energy correlates with the depletion of the first 
solvation shell and the related increase in the strength of 
the first-shell density fluctuations. Figure [6] additionally 
supports this view. Here we compare experimental [32| 
and simulated 27] atomic mean-squared displacements 
of myoglobin (small points) with our calculations of the 
diffusivity of water in the simulation box and the change 
in the number of waters in the first solvation shell. All 
parameters have been normalized to their correspond- 
ing values at 300 K to bring them to the common scale. 
The remarkable result of this comparison is that the av- 
erage number of waters in the first solvation shell fol- 
lows very closely the atomic displacements changing its 
temperature slope at the point of dynamical transition, 
Ttr = 220 K. The increased mobility of the protein is 
therefore related to the increased translational mobility 
of waters |29l l43j caused in turn by the creation of the 
high-temperature hydrophobic interface [35|. 

The diffusion coefficient of water in the simulation box 
is plotted separately vs the inverse temperature in Fig. 
[3 where we also compare our results to previous simula- 
tions by Kumar et al [6| and by Lagi et al [21] . The diffu- 
sion coefficient was calculated from the Einstein equation 
and the reported values are averaged over all waters in 
the simulation box. The different magnitudes of diffusiv- 
ity compared to previous reports @, l2l| are related to the 
different force fields used, but, more importantly, to the 
different fractions of water molecules in the simulation 
sample. Given that all molecules in the smallest sam- 
ple in Fig. [7] belong to the interface [2l| , it is not that 
surprising that these data show the slowest diffusion, in 
agreement with the common expectation of slower diffu- 
sion of waters in thin interfacial layers [44|. Neverthe- 
less, despite the use of a much larger number of waters 
{Nu, = 5886 vs N^ = 484) , we confirm here the existence 
of a crossover in the Arrhenius slope of water's diffusion 
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FIG. 7: DifFusion coefficient of TIP3P water calculated from 
all Nw = 5886 water molecules in the simulation box (circles). 
The triangles indicate the results of Ref. G for the simulation 
box containing A''^, — 1242 waters; the squares denote the 
results from Ref. [22 with N^ = 484 per simulation box con- 
taining two protein molecules. Closed and open points indi- 
cate temperatures above and below the dynamical transition 
temperature, respectively. 



coefficient observed earlier in Ref. |21L 

Two observations are relevant in respect to the dif- 
fusivity data shown in Fig. [T] First, the temperature 
law is Arrhenius both above and below the transition 
temperature with the slope decreasing at lower temper- 
atures, in accord with previous observations [alSll. Sec- 
ond, the transition temperature is shifted down to 200 
K compared to 220 K found in simulations of partially 
hydrated proteins in Ref. [2l|. The first observation im- 
plies that we observe only a change in the character of a 
secondary, Arrhenius relaxation, as indeed often seen for 
electron transfer in proteins 11], instead of a fragile-to- 
strong transition. This fact might be related to the often 
reported [45| disappearance of a relaxation in confined 
water most closely related to our simulation conditions. 
Since D{T) follows closely the decrease in the number of 
hydration-shell waters (Fig. [5]) a connection of the break 
in the slope to a secondary process produced by collec- 
tive density fluctuations of the hydration shell seems a 
reasonable explanation. The second feature might imply 
that the existence and position of the transition temper- 
ature depends on the fraction of surface waters in the 
system. While all waters in the simulation setup in Ref. 
w\ belonged to the surface, only roughly 10% of waters 
in our simulations find themselves in the first solvation 
shell (Fig. [3] and Table [l|. Likewise, we have obtained a 
fragile-to-strong crossover by considering only first-shell 
fluctuations in Fig. [51 but it is already washed out for 
the diffusivity averaged over several hydration layers. 

Our data, while pointing mostly to the interfacial ef- 
fects as the reason for the dramatic rise of A^ (T) , do not 
entirely exclude crossing the Widom line, a bulk prop- 
erty of water, from the picture. While the global rise 
of the intensity of electrostatic fluctuations within pro- 
tein is linked to the density fluctuations of the interface, 
the spike of Xw{T) at T = 220 K and the corresponding 
slowing down of the Stokes shift relaxation might well be 
linked to the crossing of the Widom line. The increased 
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FIG. 8: Reorganization energy obtained from the variance of 
the Coulomb energy gap: water component A™ [circles, Eq. 
Q], the protein component Xp from the variance of AEp 
(diamonds) and the total water/protein reorganization energy 
Atot from the variance of the total Coulomb energy gap AE^ + 
AEp (squares). 



cooperativity of water's fluctuations at this temperature 
causes a behavior similar to the critical slowing down 
with a peak in a second energy cumulant, heat capacity 
for bulk measurements |(|| and reorganization energy in 
our case. 

The spike of Xw{T), barely seen on the 1 ns obser- 
vation window, becomes more pronounced on the 10 ns 
time-scale, as is shown in Fig. [SI where the reorganiza- 
tion energies for water, protein, and the full reorganiza- 
tion energy from water/protein electrostatic fluctuations 
were collected from the entire 10 ns trajectories by cal- 
culating the variance of the total Coulomb energy gap 
AE^ + AEp. This observation suggests that some col- 
lective motions of water, signiflcantly cut off on the 1 
ns time-scale, contribute to the peak and become more 
pronounced on a longer observation scale. 



V. CONCLUSIONS 

In conclusion, we have found a dramatic increase in 
the breadth of water-induced electrostatic fluctuations 
inside the protein with increasing temperature. We link 
this increase to the creation of the hydrophobic interface 
at extended hydrophobic patches of the protein. What 
has escaped the attention of all studies of the dynami- 
cal transition in biopolymers is the onset of hydrophobic 
solvation occurring at the same temperature Ttr as the 
dynamical transition. It might be true that the creation 
of the hydrophobic interface with its large extent of den- 
sity fluctuations and intense electrostatic noise is closely 
linked to the dynamical transition, although we do not 
currently have any additional data supporting this view. 
However, if this view is correct, there should be a criti- 
cal polypeptide dimension below which the macroscopic 
hydrophobic interface does not form [35| and no dynam- 
ical transition exists. In fact, v ery recent measurements 
of terahertz dielectric response [40| of hydrated polypep- 
tides of different lengths have indicated the existence of 
an exactly such critical polymer length below which the 
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dynamical transition disappears. 

Wc found that the dynamics of electrostatic fluctu- 
ations are coupled to fast /3 relaxation of the hydra- 
tion shell. The redox activity of proteins can there- 
fore be classified as hydration-shell-coupled, according 
to the classification suggested by Fenimore et al [33 |. 
Although this coupling carries similarities with aqueous 
mixtures of simple glass- formers Q, proteins are not just 
large molecules. The formation of the hydrophobic inter- 
face is related to a particular lengthscale of hydrophobic 
patches (~ 1 nm 35]) which does not exist for small 
hydrated molecules. Not surprisingly, large-amplitude 
electrostatic fiuctuations observed here are not usually 
seen inside small molecules [l^, although this feature 
might extend to other patchy hydrophobic surfaces, such 
as lipid membranes and dendrimeric structures. 



It remains to be seen whether and how the gigantic re- 
organization energy found at high temperatures is related 
to the biological function of metalloproteins belonging to 
energetic electron-transfer chains. One can anticipate, 
from a general perspective, that a significant increase in 
the amplitude of electrostatic fluctuations can help in re- 
ducing barriers for chemical transformations by allowing 
better chances for favorable configurations from a broad 
fluctuation spectrum. 
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